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ABSTRACT 

Aims. We revised the treatment of interstellar dust in the KOSMA-r PDR (photo-dissociation region) model code to achieve a con- 
sistent description of the dust-related physics in the code. The detailed knowledge of the dust properties is then used to compute the 
dust continuum emission together with the line emission of chemical species. 

Methods. We coupled the KOSMA-r PDR code with the MCDRT (multi component dust radiative transfer) code to to solve the 
frequency-dependent radiative transfer equations and the thermal balance equation in a dusty clump under the assumption of spheri- 
cal symmetry, assuming thermal equilibrium in calculating the dust temperatures, neglecting non-equilibrium effects. We updated the 
calculation of the photoelectric heating and extended the parametrization range for the photoelectric heating toward high densities 
and UV fields. We revised the computation of the H 2 formation on grain surfaces to include the Eley-Rideal effect, thus allowing for 
high-temperature H2 formation. 

Results. We demonstrate how the different optical properties, temperatures, and heating and cooling capabilities of the grains influ- 
ence the physical and chemical structure of a model cloud. The most influential modification is the treatment of H2 formation on grain 
surfaces that allows for chemisorption. This increases the total H 2 formation significantly and the connected H 2 formation heating 
provides a profound heating contribution in the outer layers of the model clumps. The contribution of polycyclic aromatic hydrocar- 
bons (PAH) surfaces to the photoelectric heating and H 2 formation provides a boost to the temperature of outer cloud layers, which is 
clearly traced by high-/ CO lines. Increasing the fraction of small grains in the dust size distribution results in hotter gas in the outer 
cloud layers caused by more efficient heating and cooler cloud centers, which is in turn caused by the more efficient FUV extinction. 



1. Introduction 

Photo-dissociation regions (PDRs) are interstellar neighbors of 
HII regions that are already shielded from extreme ultraviolet 
(EUV) photons with an energy sufficient to ionize atomic hy- 
drogen, but where the physical and chemical conditions of the 
atomic and molecular interstellar medium (ISM) are still gov- 
erned by the r emaining far-ultraviolet (FUV ) radiation of nearby 
massive stars dHollenbach & Tielenslfl999h . Here, the transition 
from the atomic to the molecular ISM takes place, giving rise to 
a rich astrochemical environment. The combination of a chem- 
ically rich pool of species that includes a large reservoir of free 
electrons and strong energetic excitation conditions produces a 
wealth of spectroscopic emissions, carrying information on the 
internal chemical and physical structure of the PDR. 

Comparing the results of theoretical and numerical calcula- 
tions to the observed spectral line emission of PDRs is a fre- 
quently used method to infer their local physical conditions. 
Numerical models to describe the physical and chemical pro- 
cesses in molecular clouds have been successfully used for many 
years. They involve simultaneously solving the problem of a) 
radiative transfer for a cloud of gas and dust; b) the chemical 
structure by balancing formation and destruction processes for 
all included species; c) determining the thermal structure of the 
cloud by balancing all important heating and cooling processes. 

Properly taking into account all relevant chemical and phys- 
ical processes is extremely time-consuming in terms of com- 
putational power and usually leads to a large number of sim- 
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plifying assumptions in the treatment of the physics and the 
chemistry. The exact choice of these simplifications depends on 
the field of application, on the expertise of the modeller, and 
on the available physical and chemical knowledge. Many pro- 
cesses are still not fully understood if not largely unknown. A 
prominent exa mple is the description of interstellar dust (see 
lDrainell2003al for a review on interstellar dust). Many impor- 
tant physical and chemical processes require a good knowledg e 
of the properties of interstellar dust grains (lAbel et al.l l2008). 
Chemical reactions on the surface of dust grains appear to be 
the only efficient for mation route for a number of importan t 
astrochemical species (iGarrod et alJl2008l: lHall & MilrarlhoiOT) . 
The specifics of their reaction kinematics depend to a large de- 
gree on the material properties and the structure of the dust 
grains. Porous, spongy grains provide a large surface for chem- 
ical reactions but might h inder the release of the newly formed 
species into the gas phase dLeger et al.ll985l : lRoberts et al.l2007b 
iTaquet et al.l I20121) . Not only the material and shape, but also 
the size of the grains might be important. Very small grains 
and very large molecules, such as polycyclic aromatic hydro- 
carbons (PAHs), efficiently contribute to the gas heating by 
means of the photoelectric (P E) effect dBakes & Tielens|[l994t 
IWeingartner & Draindl2001al) while larger grains are dominat- 
ing the scatt ering and absor ption of FUV photons. For a review 
on PAHs see lTielensl(l2008l) . 

The rapid development in observational techniques over the 
last two decades left us with a vast amount of spectroscopic data 
that defies reproduction with simple numerical PDR models. We 
are far away from bein g able to really und erstand the conditions 
in any observed PDR dRollig et al.ll2.007h . In this paper we de- 
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scribe a number of steps that have been taken to increase the 
modeling power of the KOSMA-t PDR model code. We tried to 
evolve the code toward a consistent treatment of dust physics and 
chemistry. By calculating the optical and energetic properties of 
a model cloud with a given dust composition we are now able 
to describe the full spectral emission characteristic of a model 
cloud including dust continuum and line emissiorQ. 

2. Code description 

2. 1. The KOSMA-t PDR model code 



The KOSMA-r code dRollig et al.l|2006l) . a welltested and ma- 
ture PDR model code dRollig et al.ll2Q07h . applies spherical ge- 
ometry to the problem of simultaneously solving the chemical 
and energy balance in an interstellar molecular cloud. KOSMA- 
r is equipped with a modular chemical network, i.e., chemical 
species can easily be added or removed from the network and 
the network will rebuild dynamically, including isotopic chem- 
istry of 13 C and 18 0. After obtaining the physical and chemical 
structure of the model cloud, a radiative transfer code is applied 
to calculate the resulting emissivities for the spectral line emis- 
sion. Because of the spherical geometry of the model clump, 
radiation can reach a point inside the clump from any direction. 

2.2. Multi-component dust radiative transfer model 

The multi-component dust radiative transfer (MCDRT) code al- 
lows one to solve the frequency-dependent radiative transfer 
equations and the thermal balance equation i n a dusty clu mp 
under the assumption of spherical symmetry dYorkdfl980l) . A 
detailed description of the code's main features is given in 
ISzczerba et alJ (1 1 997b . The code includes isotropic scattering, 
which is important at high optical depths. Further modifications 
to that version of the code have been made. First, we have added 
the possibility to solve the radiative transfer simultaneously for 
a number of dust sorts (NDS) i that may or may not be co- 
spatial, and obey any dust size distribution. The grain size dis- 
tribution is given by the general relation valid for grains with 
radius _ < a < a* + ; 



drii(a) — nfi(a)da, 



(1) 



where: n,(a) is the number density of grains with size < a and n 
is the number density of H nuclei (in both atoms and molecules: 
n — riu + 2«h,)- The function has a modular structure and 
can be changed easily to a different form. Second, because we 
adapted the code to the conditions typical for PDRs it was neces- 
sary to introduce an option for switching off the central source. 
In addition, because intense sources of an interstellar radiation 
field are ubiquitous in PDRs, we changed the spatial distribu- 
tion of impact parameters (originally logarithmically spaced be- 
ginning from the inner shell of the clump) along which the ray 
equations are integrated. The end result of this code is the in- 
frared radiation flux emitted by the clump, which can be com- 
pared to the radiation observed from real PDRs. However, as a 
byproduct, the code provides the mean intensity of the radiation 
field, Jxir) and the dust temperature, T^ \{r, a), which for now is 
computed from the assumed thermal equilibrium for each dust 
species and dust grain size at each radius r of the clump. These 
two quantities (mean radiation field and dust temperature) are 



now used in the KOSMA-t PDR model code to obtain a self- 
consistent gas-dust model. 

2.3. Dust models 

For the purpose of th is paper we selected th r ee dust models of 
interstellar dust from IWeingartner & Draind d2001bl) (hereafter 
WD01), and the MRN dust model dMafhis et al.ll 1977b for com- 
parison. The MRN interstellar dust model was constructed by 
fitting a dust model to the Galactic mean extinction curve. It 
consists of two separate dust populations, one for "astronomi- 
cal silicates" (sil) and one for graphite (gra). In the MRN model 
the grain size distribution is identical for both dust components 
and is given by 



dn\(a) - nA\a 3 5 da, 0.005 < a < 0.25/j.m, 



(2) 



1 Other PDR model codes are also able to perfor m wavelength- 
depen d ent UV radiative transfer. Compare for example iLe Petit et alJ 
d2006h : lAbeletalJd2005h . 



where i stands for sil and gra, A s n = 10 25 10 , and A gra = 10 25 13 
(WD01). 

Recently gathered observational constraints led WD01 to 
construct new models for interstellar dust. The authors proposed 
that interstellar dust is composed of four main components: 
astronomical silicates, graphite, and populations of very small 
grains that consist of neutral (PAH ) and ionized (PAH + ) parti- 
cles. The grain size distributions for these dust components fol- 
low the general definition given by Eq. (Q} and their functional 
forms (fi(a) in Eq. ([TJi) are given by Eqs. (2-6) of WD01. The 
authors employ a power-law form for the large grains and a log- 
normal distribution for the very small grain population. Free pa- 
rameters of the fi(a) functions were determined by the best fit 
to the average extinction with Ry = Ay/E(B - V)) of 3.1, 4 and 
5.5 (see WD01 for details). Ay is the absolute visual extinction 
at V = 5500 A and E(B - V) = A B - Ay. The three dust mod- 
els f rom WD01 correspond to line s 7, 21, and 25 from Table 
1 in IWeingartner & Drainel d2001bh and are denoted WD01-7 
(Rv = 3.1), WD01-21 (R v = 4.0), and WD01-25 {Ry = 5.5). 

To compute the extinction efficiency Q ext and albedo u> for 
each dust compone nt, we assumed that grains are spherical and 
used the Mie theory dBohren & Huffmanll983l) . For sili c ates an d 
graphite we used the dielectric constants from iDraind d2003bb . 
while for very smal l grains we followed the approach given by 
iLi & Drainel d200ll) . The minimal and maximal grain size for 
the size distribution of each dust component are given in Tab. [4] 
For very small grains we used 17 grain sizes, while for bigger 
grains we divided each dex of grain sizes into ten sizes, keeping 
equal distance in log(a). The wavelength coverage used in the 
MCDRT code extends from 10 A to 3000 fim and is split into 
333 wavelengths. 

3. Impact on the PDR model 

In the following we describe how the new, full radiative trans- 
fer (RT) computation compares to the old approximation used 
in the KOSMA-t model. For this purpose we evaluated a ref- 
erence model clump with varying FUV radiative transfer and 
different dust models. We kept the following model parameters 
constant for all models: the total clump mass M =10 M , the 
total surface gas density n = 10 5 ctrT 3 , and the total F UV field 
X = 1000 in units of the Draine field dDraindll978h . We ap- 
plied a power law density gradient with power index 1.5 with 
a constant central gas density for radii smaller than Rt t/5. This 
implies R tot = 2.46 x 10 17 cm = 0.08pc. We assumed a total cos- 
mic ray ionization rate of molecular hydrogen £w = 5 x 10" 17 
s _1 . The chemistry is based on the UMIST 2006 database for 
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astrochemistry UDf^Q dWoodall et al.ll2007l) . using 47 chemical 
species and 490 reactions in total. 

3.1. UV continuum radiative transfer 
3.1 .1 . Dust extinction 

For any given line of sight, the optical depth at F UV wavelengths 
tfuv = ct^fuv^h ^Sternberg & Dalgarnoll989h . with <r diFUV be- 
ing the effective dust grain extinction cross section in the FUV, 
and the total hydrogen column density A^h- In the previous ver- 
sion of the KOSMA-t code, the local radiation field was cal- 
culated by computing the optical depth as a function of angle 
and then performing an angular averaging of the resulting inten- 
sity over the full solid angle. A detailed descriptio n of the FUV 
transfe r in the spherical model has been given by IStorzer et alj 
(Il996h . 

For a constant dust-to-H ratio the total extinction cross sec- 
tion per H nucleus can be derived from 



0"ext(/i) 



TextW 



NDS 

I 

!=1 



fi(a)C' ext (A, a )da, 



(3) 



where C' ex( (A, a) is the extinction cross section of a single dust 
particle with size a at wavelength A. Often the absolute extinc- 
tion in magnitudes is used 



A(A) = 2.51og( e )T ext W * 1.086T ex tW). 



(4) 



By introducing the average mass extinction coefficient for 
each dust component, i.e., the mean extinction cross section per 
dust mass, 



J ni(a)C' exl (A)da J ni(a)C' exl (A)da 



(K l t (A)) a 



J ni(a)4-/3Ka 3 pt,jda 



PdusU 



(5) 



where pb.; is the bulk dust grain density (i.e., the density of mate- 
rial of which dust of given sort is composed) and pd us t,i is volume 
density of the r'-th dust sort, the optical depth at any given point 
of the model clump can be written as 



„ NDS 

J dr£ (K' e 



tW)a Pdust.i 



Relating this to the gas density provides 



Text (A) = m H 



/NDS 
, Pdust V 
dr n > 



Pdust.i 
Pdust 



(6) 



(7) 



where n oth = p ga s, and niu is mass of a single H atom. Using 
*P = Pdust/ (h»ih) as the total dust-to-H mass ratio, and ¥ = 

NDS NDS 

2 = 2 Pdust.i /(«'«h) we can write 



with 



T e *(A)=m H J drn^¥{K ext (A)), 

, NDS 



(8) 



(9) 



total MRN 

total WDO 1-7 

- - silicates 

- - graphite 
- PAH°+PAH + 




Fig. 1. Extinction cross section per H nucleus for the WD01-7 
dust model (dashed lines) and for the MRN model (solid lines). 
The total cross sections are represented by thick lines. The con- 
tribution from each dust component is shown by red lines for 
silicates, blue lines for graphite, and green lines for very small 
grains. 

For a constant dust-to-H mass ratio the integration over r in 
Eq. (O can be performed giving an H column density A^h using 

Text(^) = A^HWxtU)), (10) 

i.e., o- ext (A) = (K^tiA^mn^. 

3.1 .2. Influence of the dust composition 

Figure Q] shows the total extinction cross section per H nu- 
cleus computed for the WDO 1-7 (thick dashed line) and the 
MRN model (thick solid line) of interstellar dust. The other 
lines show contributions from each dust component incorporated 
into the WD01-7 model (silicates - dashed, red in the electronic 
version; graphite - dashed, blue in the electronic version; and 
PAH°+PAH + — dashed green line in the electronic version), 
and the solid lines present contribution from silicates (red) and 
graphite (blue) in the MRN model. One can see that in the WD01 
model populations of very small grains (PAH and PAH + ) are re- 
sponsible for a large part of the total extinction, especially for 
the 2200 A band, while large grains dominate the total mean 
Galactic extinction in the MRN model. 

Fig. [2] shows A,\IAy for the dust models tested in this pa- 
per: MRN (solid), WDO 1-7 (long-dashed; red in the electronic 
version), WDO 1-21 (short-dashed; green in the electronic ver- 
sion), and WDO 1-25 (dot-dot-dashed; blue in the electronic 
version). The FUV range relevant for photoelectric heating 
and photo-dissociation reactions extends between 13.6 and 6 
eVQ. Consequently, we define an FUV-to-V color as &fuv = 
(A(A)/A(V))a = (A,i/Ay),i where the averaging is performed 
over an energy from 6 to 13.6 eV and Ay is the visual extinc- 
tion. The vertical dotted lines in Fig. [2] indicate this range of 
averaging and the horizontal lines show fcpuv for the considered 
models. Note that the MRN and WD01-7 dust models have the 
same value of fcpuv = 3.339. All values of the average relative 
dust extinction are collected in Tab.|4] From Fig.[2]it can be seen 
that £puv is significantly decreased in the WDO 1-21 and WD01- 
25 dust models. This is consistent with the higher values of Ry, 
i.e., a shallower gradient of A(A)/Av toward shorter wavelengths 
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H 2 only absorbs FUV photons via Lyman and Werner electronic 
transitions in the 912-1 100 A range. 
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Fig. 2. Relative extinction for the dust models of MRN (solid) 
, WDOl-7 (long-dashed; red in electronic version); WDOl-21 
(short-dashed; green in the electronic version) and WDOl-25 
(dot-dot-dashed; blue in the electronic version) of interstellar 
dust. The vertical dotted lines at wavelengths of 912 and 2066 
A, corresponding to 13.6 and 6 eV, show the range of averaging. 
The horizontal lines show fcpuv for the considered models. 




Fig. 3. Mean intensity Jx (f) for a model clump of M =10 M , 
n — 10 5 ctrT 3 , and^ = 1000. The lines show Jx(r) for radii in 
steps of 10%offl tot . 

causing a lower UV extinction for the same Ay (see Fig. [2j». This 
is the main source for the different UV attenuation properties of 
different dust types. 

3.1 .3. Influence of radiative transfer and scattering 

As a result of the new radiative transfer computations we ob- 
tain the full, wavelength-dependent FUV radiation field Jx(r) = 
^ f hd£l at all clump radii, where I a is the specific intensity 
at wavelength A. In Fig. [3] we plot Jx(r) for radii in steps of 
R lot 1 1 for the WD0 1 -7 model . The effect of the prominent 2200 
A "bump" is plainly visible in the enhanced reduction of Jx with 
decreasing radius. 

The FUV optical depth 

Tpuv = ^FUVT"V = A^hO"D,FUV (11) 

can be described by the effective dust extinction cross section 

0d,fuv = k FVV m H ¥ (K ext (V)} . (12) 

For adopted dust models, the obtained values of <x D fuv are col- 
lected in Tab. g] 



The original KOSMA-r code did not include any angular de- 
pendence of the dust scattering but implicitly assumed pure for- 
ward scattering dSternberg & Dalgarnoll 19891) in terms of an ex- 
tinction coefficient cr D FUV . As long as the models use the same 
value of <T D Fuv, the radiative transfer should always give the 
same scaling of the UV intensity with exp( -/VhO"d,fuv ) at high 
optical depths because the radiation from scattering at other an- 
gles is quickly damped because of the longer optical path, but 
the scattering will in crease the intensity close to the surface 
(lFlannervetalJll980h . 

To study the influence of FUV scattering we performed a set 
of model calculations with different scattering properties but the 
same effective value of o- D FUV = 1.76 x 10 -21 cm 2 : 

1. The old KOSMA-t FUV calculations using pure forward 
scattering. 

2. The MCDRT result for a MRN dust distribution providing 
the same ctd,fuv, but where the albedo a> is artificially set to 
zero. 

3. The MCDRT result for the MRN distribution with a full 
treatment of scattering and absorption. 

In Appendix [A] we discuss how the assumption of isotropic 
scattering compares to more realistic scattering properties and 
demonstrate that it poses a clear improvement compared to the 
pure forward scattering case used in our previous model. 

Figure [4] shows a comparison of the depth-dependant mean 
FUV intensity scaled by the unattenuated mean FUV intensity 
J(r)/xo for these three models. J(r) is the total mean intensity, 
i.e., averaged over the full solid angle and integrated over the 

full wavelength range: J(r) = C^X^ Jx( r )dA. xo is the corre- 
sponding total mean intensity of the FUV field in the absence 
of the molecular cloud, i.e., the full unattenuated photon flux 
coming from An solid angle. For r = R tot , models 1) and 2) are 
slightly larger than 0.5, which would be the theoretical value for 
a fully opaque and very large clump. Higher values indicate the 
contribution from angles slightly larger than 90°. By contrast 
model 3 shows a value of 0.58, indicating significant contribu- 
tions to J{r) by scattered photons. Consequently, the mean inten- 
sity is higher for model 3) throughout the whole clump because 
of a large amount of scattered photons. The small differences 
between models 1) and 2) result from the full wavelength treat- 
ment compared to the attenuation by an average Tpuv only. For 
high optical depths models 2) and 3) converge to the same inten- 
sities as all sidewards scattered ph otons turn insignifica nt there 
because of the longer optical path dFlannerv et al. fl980l) . 

Eq. (TT~2T > turns out to be a very good approximation when 
only forward-scaterring occurs, which is clearly reproduced by 
the detailed radiative transfer computations with the albedo set 
to zero. However, only the full radiative transfer treatment can 
produce the full spatial and spectral knowledge of the radiation 
field and the proper accounting for scattered photons. 

The influence of dust composition on the resulting depth- 
dependent FUV intensity in the clump can be seen in Fig. [5] 
For the FUV intensity in the WD01 models, the following is true 
throughout the whole clump: WD0 1 -7 < WD0 1 -2 1 < WD0 1 -25 . 



3.1.4. Influence on photo-reactions 

The local UV intensity acts on the chemistry of the clump by 
heating the gas and dust and by inducing ionization and dis- 
sociation reactions. This photo-process j with a wavelength- 
dependent cross section o-j(A) proceeds at a rate dRoberge et al.l 
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Fig. 4. Mean intensity of the UV radiation vs. cloud depth. The 
black line corresponds to the old KOSMA-t result for pure ab- 
sorption. The blue line corresponds to the MCDRT result for an 
MRN distribution with full treatment of scattering and absorp- 
tion. The dashed line is the result for an MRN distribution with 
the albedo artificially set to zero. 



Fig. 6. Fit results for yf fjl for HCO (black), CH 2 (red), 
CN (green), and N, (blue ) for all dust models from 
IWeingartner & Draind (l200"lbl their Tab. 1). The abscissa gives 
the number of the dust model. The large points show the original 
values of jj for the respective molecules. 
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Fig. 5. Mean intensity of the UV radiation J(r) scaled to the un- 
shielded (empty space) Draine field xo vs. cloud depth. Different 
lines denote different dust compositions. 



Il991h 



Fj(r) = An 



J,i H 



(A)dA. 



(13) 



7i(r) is the mean intensity of radiation at radius r in photons 
cm -2 s _1 nm"'. The integral runs from Ah — 912 A to the thresh- 
old wavelength Aj for the process J3- Astrochemical databases 
often provide parametrized reaction rate coefficients to allow the 
calculation of these photo-reaction rates. For instance, UDfA 
gives the rate coefficients in the form 



r/(A v ) =xoajexp(-yjAv) 



(14) 



assuming an attenuation of the radiation provided by a standard 
MRN dust model in a plane -parallel configuration using the ra- 
diation transfer results from lFlannerv et al. ( 1980) parametrized 
in terms of the perpendicular visual extinction Ay. xo is the FUV 



4 A database of photo-ionization and photo-dissociation cross 
sections for common astrophysical species can be found at 
http://home.strw.leidenuni v.nl/~ewine/phot o/ For more details see 
van Dishocck etal. (2006) and van Hemert & van Dishoeckl ( p2008h 



Fig. 7. Comparison of explicitly fitted for HCO (black), 
CFT (red), CN (green) , and N (blue) for all dust models from 
IWeingartner & Draind d2001bl their Tab. 1) with scaling rela- 
tionship Eq. [16] The large points show the original values of jj 
from UDfA, the open diamonds show the values calculated using 
an MRN dust model, and the open squares denote the explicitly 
calculated values assuming 1 /R 



760A 



3.62. 



field strength at the edge of the cloud and the unshielded, free- 
space rate coefficient <xj assumes a mean FUV intensity of unity 
(in units of the Draine field) and FUV photons coming from all 
directions. At the edge of an optically thick molecular cloud the 
rate coefficient is about half of this value, depending on the dust 
scattering properties. 

The unshielded rate coefficient aj does not depend on the 
dust properties and we concentrated on the attenuation proper- 
ties, parametrized in form of the . Changing either the dust 
properties, i.e., the wavelength-dependant FUV attenuation, or 
the spectral distribution of the FUV field leads to a different 
result for Fj and consequently to a different fit parameter jj. 
However, performing the wavelength-dependent integration of 
Eq.[l3]for every reaction and every radial point is often beyond 
the scope of PDR models, e.g., because it is too time-consuming 
or because neither dust properties nor ionization cross sections 
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<Tj(A) are known accurately enough to justify the computational 
effort of a full integration of the wavelength-dependent radiative 
transfer. Furthermore, knowledge of the wavelength-dependent 
o-j is not available for all astrophysically relevant species and 
models have to apply Eq. [14] even if they can perform the full 
integration in principle. 

Here we derive a scaling relation between the yj values pro- 
vided by a data base like UDfA and corresponding fit values for 
a particular dust model yj — > y® to calculate the rate coefficient 

(Ay) for species j and dust sort D. For a direct comparison 
with UDfA values we started from a model cloud with a mass of 
10 3 M and a radius of 1.71 pc, large enough to neglect the ef- 
fects of the spherical symmetry close to the surface, isotropically 
illuminated and calculated the photo-diss ociation and photo- 
ioniza ti on rates of 72 species provided by Ivan Dishoeck et aD 
d2006l) : lvan Hemert & van Dishoeckl (120081) for all 25 dust mod- 
els D from Weingartner & Draine ( 2001 rj) according to Eq. Qj] 

To this end, we performed least-squares fits to 
rj(A v )/Fj(0) - exp(-y® fit Ay) up to an Ay = 10 to deter- 
mine new y D jf it coefficients for all D and j. Figure [6] shows the 
new y D jj it for four example molecules, HCO, CH 2 , CN, and N 2 
for all 25 dust models. 

The big dots on the left show the original value yj from 

UDfA0 

The figure shows that the new y , .. show the same be- 
havior along the abscissa for all four example species. There 
are also three branches of y^„, visible as prominent quantita- 
tive steps in Fig.|6]that correspond to the respective Ry values. 

The next step is to identify what dust model property is best 
correlated with J^j it - An extensive analysis of the dust properties 

reveals that we can find a wavelength A = 760A for which 



1/Ra = 



A(V) 



A(A)-A(V) 



(15) 



shows a monotonic ordering for the J^j it of the 25 dust models. 
Consequently, we chose 1 /Rj 60 ^ as the parameter describing the 
dust. TableQ]gives 1 /R D . for all 25 WD01 dust models. Using 

1—10 ' 760A & 

this parameter, we can approximate the dependence of y® on yj 
and the dust model by 



W\ + W2 



760A > 



+ W3 



K R D . 

v 760A > 



+ W4 



,R D . 

760A > 



(16) 



3 4 VJ 

+w 5 yj + wey } + wry } + ^8^75- 



Wg 



760A 



R D ' 

760A 



Equation[T6]reproduces the explicitly calculated attenuation 
of the rate coefficients within an accuracy of <x = 0.05, and a 
maximum absolute residual of 0.22. This is comparable to the 
uncertainties of J?j it due to fitting T 7 (Ay) /T 7 (0) up to a differ- 
ent maximum Ay. Note that for the fitting we removed seven 
species (photo-dissociation of CH + , SH + , OH + , HCO + , CO, 2 , 
and SiO ) from the data set because their y® represent strong 
outliers with respect to the remainder of the data set. We discuss 
possible reasons for these deviations in Appendix iBl To apply 
Eq. [16] to them, we had to use corrected coefficients fj instead 
of the UDfA vales when computing y?. The fj values are given 
in Tab. [3] 



5 Ivan Dishoeck et alJ i2006h and lvan Hemert & van Dishoeckl d2008h 
provide their own fitted values of y^EvD (up to an maximum A v = 3) for 
all species for which they also give tables of <x,(/l). 



Tabl e 1. Dust model parametrizati on of all WD01 dust models. 
See lWeingartner & Drainel(l2001bl) for more details. 
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25 
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Notes. 

(n) Dust model D, see Tab. 1 in WD01 
(6) Ry = A(V)/£b-v> rauo of visual extinction to reddening. 
(c) C abundance in double log-normal very small grain population. 
(<i) Case A: variable grain volumes. Case B: grain volumes fixed at ap- 
proximately the values found for R\i = 3.1. 
^ l/R A = (A(V)/(A(A)-A(V)))- [ 

{ - f) Ordering of the various dust models from lowest to highest value of 

1/^760A- 



Table 2. Fit parameters in Eq. [16] 



fit parameter 



W\ 


1.39460 


"'2 


-0.27655 


W3 


0.02053 


>1'4 


-0.00295 


"'5 


-0.45990 


We 


0.08725 




-0.01616 


Wg 


0.24747 


Wt) 


-0.01677 



It turns out that the original 

? MRN 



MRN distribution with 



l/R MK 2 = 4.03 does not exactly follow the parametrization 



of Eq. \ T6\ based on the WD01 dust models, but that it can be 
easily used in the same function when adjusting the parameter 
1//?mrn = 3_ 62 . Figure [7] compares yf fu for HCO, CH 2 , CN, 
and N 2 with the results of Eq.[16] In addition we show the MRN 
results in the plot: the original yj, the explicitly calculated y^ 1RN , 
and the -y^ 1 ™ assuming 1/^^^ = 3.62 using large points, open 
diamonds, and open squares, respectively (the arrows in the fig- 
ure demonstrate the shifting of the rescaled y® closer to Eq.fTSTl. 
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Table 3. Corrected fj for the outlying species. 



species y ; yj 



CH + 


2.5 


2.11 


SH + 


1.8 


1.39 


OH + 


1.8 


2.76 


HCO + 


2. 


3.036 


CO 


2.5 


3.03 


o: 


2. 


1.61 


SiO 


2.3 


1.87 
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Fig. 8. Dust temperatures of silicate grains for a model clump 
of M =10 M Q , n = 10 s cm" 3 , and* = 1000 using the WD01-7 
dust model. The different lines denote different grain sizes. 

It is unclear why the jj from UDfA, i.e., the MRN values, 
do not follow the parametrization of Eq. fTo] There i s a num ber 
of possible reasons for the deviation. IWoodall et a D (f2007h did 
not use a Draine FUV radiation field but assumed a 10000 K 
black-body radiation field. The different spectral shape usually 
accounts for approximately 10% difference but can in a few spe- 
cial cases result in a significant difference. Secondly, different 
^v.max up to which the jj is fitted give different values of y ; ( see 
Appendix iBl for a more detailed discussion). The details of the 
radiation scattering is another possible reason for different re- 
sults. Assuming different anisotropy factors g and dust albedo co 
will change the radiation field and accordingly affect the jj (see 
Appendix lAl for a discussion on the effect of different values of 

We conclude that the scaling relation Eq.[T6]is able to rescale 
the tabulated jj from UDfA to the dust-specific y?„. ( (fitted to the 
fully calculated r j(Ay) /T ;(0)) with an accuracy of about 10%. 

3.2. Dust temperatures 

Dust temperatures are calculated independently for each dust 
component and size. Figures [8] [9] and [TOl show depth-dependent 
grain temperatures for a number of differently sized grains for 
three of the four components. The plot for the ionized PAHs is 
omitted because their temperature behaves like that of neutral 
PAHs. The general behavior is that smaller grains exhibit higher 
grain temperatures. 

Figure QT| compares grain temperatures of the components 
of the MRN and WD01-07 dust model at the surface of a model 
clump as a function of grain radius a,. Graphites tend to have 
a higher temperature than silicate grains except for very large 



the KOSMA-r PDR code 
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Fig. 9. Dust temperatures of carbonaceous grains for a model 
clump of M =10 Mo, n = 10 5 cm 4 , and* = 1000 using the 
WD01-7 dust model. The different lines denote different grain 
sizes. 



100 




0.001 0.002 0.005 0.010 0.020 0.050 

cloud depth [pc] 



Fig. 10. Dust temperatures of neutral PAHs for a model clump 
of M =10 M , n — 10 5 cm -3 , and* = 1000 using the WD01-7 
dust model. The different lines denote different PAH sizes. 



grain sizes. Small grains tend to be hotter than bigger grains. 
The PAHs are heated much more efficiently than the larger grains 
and have significantly higher grain temperatures. For simplicity, 
we deliberately neglected the non-equilibrium heating of very 
small grains and PAHs here and assumed equilibrium heating in 
this paper. In a subsequent paper, we will include the compar- 
ison with real observations accounting for the stochastic heat- 
ing of the very small particles. The dependence on a,- vanishes 
for sufficiently high extinctions i.e., once most of the remaining 
photons have wavelengths longer than the grains. The resulting 
equilibrium grain temperature then only depends on the absorp- 
tion coefficient of the material. 

The detailed knowledge of the dust temperature as a func- 
tion of radius, grain size, and grain type enables us to treat all 
dust-related processes separately for all grain sizes and types. 
However, since most relevant processes are grain surface reac- 
tions, e.g., formation of the H2 or photoelectric heating, it is 
usually sufficient to calculate surface-averaged quantities. For 
the z'-th dust component we can define the relevant mean dust 
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- MRN 




Td, silicate 


MRN 




Td, graphite 


WDOl 


-7 


Td, silicate 


WDOl 


-7 


Td, graphite 


— WDOl 


-7 


^/,PAH° 


— WDOl 


-7 


Td, PAH + 





Fig. 11. Grain temperature for the components of the MRN and 
WDOl -07 dust models as a function of the radius a, in nm. The 
grain temperatures are calculated at the surface of a model clump 
of M =10 M Q , n = 10 5 enr 3 , and^ = 1000. The different lines 
denote different dust components. Dashed lines show the MRN 
components and solid lines show the WD01-07 components. 



temperature as 



(T d ,i) a 



Jg' + Td,i(a)rii(a)a 2 da 



J^"' + rii(a)a 2 da 
Averaging over all dust sorts is done according to 



(T d ) 



NDS 



2 ( T d.i)a 



Q~dj 



(17) 



(18) 



with the geometrical cross section of the i-th sort (Td ,; and the 
total dust cross section cr rf . We will use the short-hand notation 
Td = {Tdhos and T dJ = {T d j) a from here on. 

Figure [12] shows the mean dust temperature T d j for the WD- 
7 dust model as a function of cloud depth for all four dust 
components. As at the surface (Fig. fTTT i. silicates, indicated by 
the solid blue line, are cooler than graphite grains (solid, red 
line). Graphite and silicate grain temperatures are very close, but 
the mean silicate temperature drops slightly faster for increased 
cloud depths than the graphites. The PAHs have a much higher 
Tdj. Neutral PAHs tend to be marginally hotter than the ionized 
PAHs. The black line denotes the average dust temperature Td 
for the WDOl -7 dust. 

The influence of the different dust models can be seen from 
Fig. [13] The black line shows the old dust temperatur e calcu- 
lation in KOSMA-r following IHollenbach etaTT (1 1 99 lh . which 
gives a central dust temperatures T d ~ 22 K. Accounting for the 
full wavelength dependence and the detailed grain size distribu- 
tion gives a significantly lower Td ~ 13 K at the center of the 
model clump. On the other hand, it also leads to warmer dus t at 
the surface of the cloud compared to IHollenbach et all (1 19911) . 

The dust temperature for the WDOl models follow their cor- 
responding FUV intensities (see Fig. [5]). The MRN model has 
a much lower Td- Figure QT] shows that the Tdj for MRN and 
WDOl -7 are very similar for equal grain sizes. However, the 
temperature of the PAHs in the WDOl -7 model contributes sig- 
nificantly to the average dust temperatures in Eq. ( [T81 and leads 
to higher Td in models with PAH content. 

The lower central dust temperatures of models with WDOl 
dust types also influences the gas temperatures. Models with 



Fig. 12. Mean dust temperature for the WDOl -7 dust model for 
a model clump of M =10 M , n = 10 5 enr 3 , and^ = 1000. The 
different lines denote the average temperature of the individual 
dust components. The black line shows the overall mean dust 
temperature. 





old FUV 




- - - WDOl - 7 - 




--- WDOl - 21 - 
-- WDOl -25 " 




MRN 











Fig. 13. Mean dust temperature as a function of A v for a model 
clump of M =10 M , n = 10 5 enr 3 , and X = 1000. The different 
lines denote different dust compositions. 



WDOl dust show significantly lower central gas temperatures 
because gas-grain collisions are becoming the dominant cooling 
process. 

3.3. H 2 formation 

Molecular hydrogen is the most abundant molecule in the uni- 
verse. The formation of H2 most effectively takes place on the 
surface of dust grains. The proper treatment of the H? formation 
is cru cial in any numerical PDR model (see e.g. lLe Bourlot et all 
120121 for a recent update of the Meudon PDR code.). We up- 
dated the KOSMA-r calculation of the H2 formation efficiency 
Rn 2 from the method described by SD89 to account for ph- 
ysisorption and ch emisorption binding ene r gies on silic ate and 
graphite surfaces (ICazaux & Tielensl (2002L |2004 l20ld (here- 
after CT)). 

This model is based on two main points: (1) H atoms can 
stick to the grain surface in two different binding sites, a ph- 
ysisorption site and a chemisorption site. Physisorption is a 
weak atom-surface binding resulting from a Van der Waals 
force (dipole-dipole interaction). Chemisorption is a strong 
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Table 4. Parameter of the WD01 and MRN dust models 
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atom-surface binding interaction involving the valence electrons, 
known as covalent bond. These interactions determine whether 
atoms move on the surface, and thus, how quickly they can meet 
a recombination partner. (2) The mobility of a surface H atom re- 
sults from the combination of two physical processes: tunneling 
and thermal diffusion; tunneling dominates at the lowest temper- 
atures while thermal diffusion is most important at the highest 
temperatures. The calculation o f the transmission coefficients, t o 
go from site to site, are given in lCazaux &Tielensl (12004 l2010h . 

At low temperature (T < 100 K), H2 formation involves 
the migration of physisorbed H atoms. At higher temperatures 
(T > 100 K), H2 formation results from chemisorbed H recom- 
bination. The presence of these two types of binding sites allows 
Information to proceed relatively efficiently at low and elevated 
temperatures. 

For the dust component i, the formation rate can be written 

as 

1 



^H 2 ./ = r 7n n v u n d jCr di e ni ,iS H , 



(19) 



where «h and vh are the number density and thermal veloc- 
ity of H atoms in the gas phase, n dd (Td,i is the total cross sec- 
tion of the grain compo nent i, eu 2 j is the formation efficiency 
(ICazaux & Tiel ens 2004) on surfaces of the i-th component, and 
S h is the sticking coefficient of the H atoms 



1 + 0.04 



T g + T d 
100 



1/2 



+ 0.2- 



100 



+ 0.08 



100 



2\ 



(20) 

which depends on the gas and dust temperature, T g and T d 
jHollenbach & McKed 1 1979b . The thermal velocity is 1.45 x 
10 5 (r ? /100)'^ 2 cm s _1 . The z'-th dust cross section per H nucleus 
is given by 

f° i+ 2 
o~dj = I n,(fl)7Tfl da . (21) 

kJ a. _ 




T d [K] 



Fig. 14. H2 formation efficiency for carbon and silicate surfaces. 
A and S denote the temperature regimes of the corresponding 
terms in Eq. d22l . The das h-dotted line denotes the stand ard for- 
mation efficiency given bv lHollenbach & McKed (1 1979b . 



Table Bl lists cr d for the different dust models. ICazaux & Tielensl 
(l2002h expressed the general expression for the z'-th H2 formation 
efficiency: 

en 2 ,t = - 7 - ■ (22) 
SI + 1 + s 

eH 2 ,/ is limited by three terms: the first term JA. prohibits the evap- 
oration of the newly formed molecules at low temperatures. The 
second term, equal to unity, dominates at higher grain tempera- 
tures; all incoming H atoms leave the surface as H2. The third 
term H governs the high-temperature regime, where evapora- 
tion of physisorbed atoms removes H from the surface before 
it can recombine to H2. It competes with the hopping of atoms 
into chemisorbed binding sites where the atoms can remain on 
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r e =50K 



10" 




T„ [K] 



Fig. 15. H2 formation rate Rn n as a function of dust tempera- 
ture Td for different formation models and dust compositions. 
The gas temperature T g is set to 50 K. The dashed lines in- 
dicate the-low temperature formation in the original CT for- 
mation rates with + 0. SD89 is the formation rate from 
ISternberg & Dalgarnol dl989h . 



the surface to form H2 at significantly higher temperatures. The 
numerator describing the chemisorption, terminates the H2 
formation at much higher temperatures, above several hundred 
Kelvin, when even chemisorbed atoms start to evaporate. 

Figure [14] shows en 2 as a function of Tj for silicate (blue 
line) and gr aphite (red line) surfaces. Th e old formation effi- 
ciency from IHollenbach & McKed (fl979l) is shown as a dash- 
dotted line for comparison. When applying this formation in the 
framework of the PDR model, see Fig. [T4] the term JA. turned 
out to be problematic because it shuts down the H2 release into 
the gas phase once T < 10K, accumulating all molecules on the 
grain surfaces. We question this term because Eq. (1221 does not 
account for the H2 binding energy of 4.48 eV that upon forma- 
tion is available to allow release into the gas phase even at very 
low temperatures. As a remedy we set J? = to ensure efficient 
H2 formation even on cold dust surfaces (Cazaux (2011), priv. 
comm. More details on how to calculate the H2 formation ef- 
ficiency are given in Appendix [C] The total H2 formation rate 
can then be written as 



NDS NDS 



(23) 



Fig. [15] compares the new and old H2 formation rates. The 
formation efficiency given b y SD89 follows the treatment of 
IHollenbach & McKed (1 19791) . which adopted a critical dust tem- 
perature of T cr x 65 K above which the Fb formation effi- 
ciency drops below 0.5. H2 formation becomes inefficient above 
» 100 K. The H2 formation in the CT dust model remains ef- 
ficient up to 1000 K. The model by Cazaux&Tielens allows 
a much more efficient H2 formation at higher temperatures. 
Additionally, the significantly larger dust surface in the MRN 
and WD01-7 model lead to a much higher FF formation rate 
than SD89. 

While accounting for chemisorbed H atoms significantly in- 
creases the H2 formation efficiency at higher surface tempera- 
tures, it is the total available grain cross section that distinguishes 
between different dust models. Table [4] also compares the total 



WD01-07 
WD01-21 
WD01-25 



10-18 r r s =50K 



10" 




Fig. 16. H2 formation rate as a function of dust tempera- 
ture Tj for different dust compositions. The gas temperature T g 
is set to 50 K. The dashed lines show the formation rate with 
suppressed FF formation on PAHs. The solid lines show the for- 
mation rate allowing for FF formation on PAH surfaces. 



WD01-7 model 



SD89 

--- C&T(2002) 

C&T(2002) no PAH's 



6 All model results accounting for CT H 2 formation shown in this 
paper assume J4 = 0. 



Fig. 17. H2 formation heating rate for a model clump of M =10 
M Q , n = 10 5 cm -3 , and x — 1000. The different lines show 
the effect of different H2 formation treatment. All models have 
been computed using the WD01-7 dust properties regarding the 
radiative transfer and the PE heating. 



dust cross section per H nucleus for the different dust compo- 
sitions. The PAH cross sections contribute the majority to the 
final H2 formation rate. The resulting H2 formation rate for one 
H atom cirT 3 is shown in Fig. [16] Solid lines show the formation 
rate as a function of T c / for a constant T g — 50 K for the WD01-7, 
WD01-21, and WD01-25 models allowing for H2 formation on 
the surface of PAHs. The dashed lines show the formation rates 
if the PAH surfaces are excluded from the H2 formation calcu- 
lation. Consequently, the dust model with the largest effective 
cross section, WD01-7 shows the highest formation rate in both 
cases. 

The strong increase in the total H2 formation rate has an 
important effect in addition to the easier H2 production. Each 
formed H2 molecule relea ses its binding ene rgy of 4.48 eV. A 
common assumption (e.g jHabart et alj|2004l) is that this energy 
is equally used to a) release the molecule into the gas phase, b) 
increase the kinetic energy of the molecule, and for c) internal 
rotational-vibrational excitation. If we assume a uniform distri- 
bution, 1 .5 eV are available to heat the gas. Here we can measure 
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the Ht formation rate across the cloud by equivalently looking 
at the H2 formation heating. 

By including H2 formation according to Eq. ( l23l l we also 
increase the formation rate significantly at low values of Ay. This 
is shown in Fig. [17] All models in the plot were computed using 
the WD01-7 dust properties regarding the radiative transfer and 
the PE he ating. The solid line shows the results for the formation 
rate from ISternberg & Dalgarnd (Il989l) . The result for the full 
H2 formation treatment from CT (with !A. = 0), suppressing and 
allowing H2 formation on PAH surfaces, are shown by the dotted 
and dashed line. The much larger surface area of the WD01-7 
dust model leads to a significant increase of the formation rate. 
Particularly the large PAH surface in that dust model contributes 
dominantly to the total grain surface and consequently to the 
total H2 formation heating. 

It is important to note that the H2 formation heating is exten- 
sive at cloud depths where most of the hydrogen is in the form 
of H. This leads to the curious situation of a strong H2 formation 
heating in the absence of H2 molecules. The gas temperatures at 
these depths can reach more than thousand Kelvin. Under these 
conditions two more H2 destruction reaction become important: 



H 2 +H — >H + H + H 



H 2 + H 2 



H 2 + H + H. 



(CI) 
(C2) 



Each destruction process requires 4.48 eV to break the binding, 
in contrast to the release of 4.48 eV during the formation. The 
much more efficient formation of molecular hydrogen by Cazaux 
& Tielens makes it necessary to adopt an additional cooling term 
in the chemical network, accounting for th e 4.48 eV energy con - 
sumption during the kinetic dissociation dLepp & Shulllll983l) . 
The cooling rate can be written as 

A kd = 7.2 x KT 12 n H: (fc|cTl»H + fcjC2l w H 2 ) erg s^cnT 3 , (24) 

where and ^ |C2l ar e tne temperature-dependent rate coef- 
ficients for reactions ( IC II ) and iC 2i . 

The role of interstellar dust in ISM physics and chemistry 
remains one of the key problems in modern astrophysics. In this 
section we demonstrate the strong effect of H2 formation on PAH 
surfaces on the chemistry and the thermal conditions in PDRs. 
However, it remains unknown whether PAHs really contribute 
to the formation of molecular hydrogen. Naturally, this is con- 
nected to the uncertainty about the detailed form and distribution 
of interstellar PAHs. 

By assuming that H formation can take place on PAH sur- 
faces the formation rate from Cazaux & Tielens predicts rates 
of a few times 10" 16 s" 1 when applied to the dust distribu- 
tions presented by WD01. This is about a factor 10 higher 
than what has commonly been found for the diffuse medium 
dJura [ll974t iBrowning et all 120031: iGrvet all 120021; IWeltv et alj 



2003 



Worfireetal]l2008h . However, the usual values of R 



2 - 4 x 10 cm s are mean values along the line of sight 
while our result s are local values and a d irect comparison is dif- 
ficult. Recently, iLe Bourlot et alj (1201 2l) presented an update to 
the H 2 formation formalism in the Meudon PDR code and also 
found a strong increase of the local H n formation rates. If we 
prevent H 2 formation on PAH surfaces, the total grain surface 
available to form molecular hydrogen is significantly reduced 
and we predict formation rates of the same order as found in the 
diffuse ISM. However, PDR models assuming diffuse gas forma- 
tion rates a re not able to reprod uce the observed extent of H, line 
excitation ( lHabart et alj 201 lT) . H, fo rmation on PAHs has also 
been proposed by Habart et aD (l2003h to explain the H 2 emission 




Fig. 18. Photoelectric heating rate for a model clump of M = 10 
M Q , n — 10 s cm" 3 , and^f = 1000. The different lines denote 
different dust compositions. Models with an M RN dust composi- 
tion ap ply the PE heating efficiency as given bv lBakes & TielenH 
(Il994l) . Dotted lines indicate a shift from heating to cooling be- 
cause of A gr > r pc . 



in p Oph, and recently iMennella et al.l (1201 2l) concluded from 
experimental studies that hydrogenated neutral polycyclic aro- 
matic hydrocarbon molecules act as catalysts for the formation 
of molecular hydrogen. 

3.4. Photoelectric heating 



IWeingartner & Draind d2001al) presented the photoelectric heat- 
ing rate for the WD01 grain size distributions in a parametrized 
form as a function of iff — (G yT)/n e , where T is the gas temper- 
ature in Kelvin and iff is in units of K'^cm 3 . To convert the FUV 
field from Habing G to Drainer units use G = 1.7 1%. These 
parametrizations are fairly accurate when 10 3 K < T < 10 4 K and 
10 2 K I/2 cm 3 < iff < 10 6 K^ 2 cm 3 . However, the parameter range 
in PDR models allows iff to r each values a s high as 10 9 K'^ 2 cm 3 
and as low as 10 _2 K^ 2 cm 3 . IWeingartnerl (120091) provided up- 
dated calculations of T pc and A ar for the dust distributions 
WD01-7, WD01-21, and WD01-25 up to iff = 10 9 K I/2 cm 3 . 
These rates are fairly well reproduced by the new analytic fits: 



C + dT c < 
pc ~ H \l + C 2 ^[l+C 3 ^] 



+ C 7 | 10" 26 ergs" 1 cm" 



and 



A gr = n e n H (D 1 (D 2 T D ' +D 4 )T D " 
( D5e D 6( ,-io.3>(-^) + j)-* 

+T Ds +D 9 T D, ° +D n ) 



(25) 
(26) 



10" 28 ergs"' cm" 3 . 

Numerical values for C, and D, are tabulated in the appendix in 
TablelD~T1andlD~2l 

Figure [18] shows the PE heating rates for a model clump of 
M =10 M , n = 10 5 cm" 3 , and x = 1000. The solid lines show 
th e photoelectric heating (P EH) rates using the method described 
bv lBakes & Tielens! (1 1994b . The PEH rate obtained for an MRN 
dust composition when solving the full FUV radiative transfer 
problem is enhanced compared to the old KOSMA-t values be- 
cause of the stronger FUV intensities due to the scattering. The 
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dashed lines in Fig. [T8l show the heating rates using the updated 
parametrizations from Eqs. ( 1251 and d26T > for the different dust 
models. Dotted lines indicate negative values, i.e., a transition 
from heating to cooling because A gr > r pc . For Ay > 1 the PEH 
rate of the WD01 dust models follow the FUV intensity of the 
respective models, i.e., the model with the highest mean FUV in- 
tensity, WD01-25, shows the strongest PEH rate while the model 
WD01-07 drops quickest. However, it is interesting to note that 
for Ay < 0.1 the PEH rate for the WD01-07 dust produces the 
strongest PEH rate, leading to the highest gas temperature of all 
models (see Fig. |T9b. 

4. Application of the single clump model 

4. 1. Impact on gas temperature 

The changes in computating the PE heating and the H2 for- 
mation, which depend on the applied dust model, have a pro- 
found effect on the thermal structure of the model clump. Figure 
[191 shows the separate effects of the changed FUV and result- 
ing PE heating (left panel) and different H2 formation treatment 
(right panel). Among the dif ferent dust models presented by 
IWeingartner & Drainel (12001 bl) the WD01-7 dust models shows 
the highest total PE heating rate because it contains the larges t 
population of very small grains (IWeingartner &Draine 2001a). 
Depen ding on iff, the PE heating rate from iBakes & Tielensl 
(fl99l can become higher than in the WD01 models despite its 
smaller population of very small grains. This is because they 
assume a higher electron sticking coefficient. In the chemically 
very active zone of 0.1 < Ay < 2, the new FUV treatment leads 
to systematically warmer gas temperatures. The net effect of the 
dust model on the PE heating strongly depends on the applied 
dust model, with the PE heating efficiency proportional to the 
size of the very small grain population. A second effect comes 
from the H2 formation heating. The right panel in Fig.[T9lshows 
how the increase in the total available grain surface raises the gas 
temperature. However, this mostly affects the outer layers of a 
model clump where the high grain temperatures require account- 
ing for chemisorbed H atoms to allow the formation of molecular 
hydrogen. 

Overall, we can distinguish two regions: For optical depths 
Ay ^ 0.5, the different attenuation of the FUV field by differ- 
ent dust distributions dominates the gas temperature. The model 
with the lowest reddening, i.e., the highest FUV intensity at 
those depths, shows the highest gas temperature. At lower op- 
tical depths, the FUV intensity is in all cases high enough that 
the total dust surface available for PE heating and H2 formation 
determines the gas temperature. Here, the models with PAHs and 
many small grains, which also produce a higher reddening, show 
the highest gas temperature. 

4.2. Impact on gas chemistry 

The described changes to the model have a significant effect on 
the chemical structure of the cloud. The influence of the assumed 
dust composition and dust size distribution affects the FUV ra- 
diative transfer and the improved treatment of dust scattering in- 
creases the FUV intensity throughout the clump compared to the 
original approximation. Chemical species that are dominantly 
formed or destroyed via photo-ionization or photo-dissociation 
processes are affected. Secondly, the stronger H2 formation effi- 
ciency leads to an increased H2 formation heating contribution, 
which in turn produces an increase in gas temperature at low Ay. 
Consequently, the abundance of the species that are dominantly 
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Fig. 20. H and H2 abundances for a model clump of M =10 
M Q , n = 10 5 cm 3 , and \ - 1000. The different lines de- 
note different treatments of the H2 formation. The solid line 
sh ows a model using the standa rd H2 formation rate as given 
by ISternberg & Dalgarnol d 19891) . The dashed and dotted lines 
show model results for H2 formation according to CT where the 
formation on PAH surfaces is switched on and off, respectively. 
All three models assume a WD01-7 dust distribution. The green 
dash-dotted line assumes an MRN dust model and uses the CT 
H2 formation. 



formed in the outer regions of the molecular clump is changed 
because the higher gas temperature accelerates their respective 
formation and destruction reactions in the gas phase. 

Both effects can go in opposite directions, but usually one 
effect dominates. In Fig. ( 1201 we show the dependence of the H- 
H2 transition zone on the applied H2 formation rate. We compare 
three models assuming a WD01-7 dust model with full radiative 
transfer and photoelectric heating treatment. The only difference 
between the model calculations is the treatment of H2 forma- 
tion. We apply (1) the stand ard H2 formation rate as given by 
ISternberg & Dalgarnol dl989l) (black, solid line), (2) the H 2 for- 
mation given by CT with H2 formation taking place on the PAH 
surface (blue, dashed line), (3) and the H2 formation given by 
CT with H2 formation on PAHs suppressed (red, dotted line). H2 
chemistry is unique in the sense that it only weakly depends on 
the gas temperature (see Eqs. ( [1911201 ). Figure[20]shows how the 
different treatment of the H2 formation influences the details of 
the H-H2 transition. In that particular model clump, the dust tem- 
perature is about 50 K and almost identical in all models, lead- 
ing to a H2 formation efficiency of unity across the dust mod- 
els shown. The only significant difference is the total dust sur- 
face of the various models. The SD89 approximation, case (1), 
has the smallest available surface (ov = 4.14 x 10~ 22 cm 2 ) fol- 
lowed by CT without PAHs and CT including PAHs (see Table 
3). This agrees with the H-H2 transition moving outwards with 
growing dust surface in Fig. [20] To emphasize this behavior, we 
added a fourth dust model to the plot, showing the results for 
the MRN model. The total surface in the MRN model is only 
slightly smaller than the WD0 1 -7 model without PAHs and con- 
sequently, both models show a very similar behavior. Please note 
that the contribution of graphite and silicate to the total dust sur- 
face is still very different in the different dust models. The MRN 
model possesses about equal surfaces in both kinds, while the 
majority of the big grain surfaces in the WD01-7 model comes 
from silicates. 
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Fig. 19. Gas temperature T g for a model clump of M =10 M , n — 10 5 cm 3 , and^ = 1000. Left panel: The different lines denote 
different dust models that affect the FUV radiative transfer and the PE heating. All models in the left panel werecalculated using the 
SD89 H2 formation. Right panel: The different lines show the effect of different H2 formation treatment. All models in the right 
panel were computed using the WD01-7 dust properties regarding the radiative transfer and the PE heating. 
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Fig. 21. CO and C + abundances for a model clump of M =10 M Q , n = 10 5 cm -3 , and^f = 1000. Left panel: The different lines 
denote different dust models that affect the FUV radiative transfer and the PE heating. All models in the left panel have been 
calculated using the SD89 H2 formation. Right panel: The different lines show the effect of different Hb formation treatment. All 
models in the right panel have been computed using the WD01-7 dust properties regarding the radiative transfer and the PE heating. 



In contrast to H2 , all other chemical species in the gas phase 
are affected by variations of the local FUV intensity and the lo- 
cal gas temperature. This is apparent, for example, in the C + to 
CO transition. Figure[2T|shows the density profile of C + and CO 
for different model calculations. In the left panel we show the 
influence of the different dust models in terms of radiative trans- 
fer and PE heating, while the right panel demonstrates how the 
different H2 formation treatment influences the chemistry. All 
models in the left panel were calculated using the SD89 H2 for- 
mation. The transition from C + to CO occurs at different Ay for 
all models. Consistent with the FUV intensities shown in Figs. 
[4] and [5] the models with the lowest FUV intensities show a C + 
to CO transition at lowest Ay, while models with a weak dust 
attenuation, e.g. WD01-25, exhibit the same transition deeper in 
the cloud. Additionally, we note an enhanced CO density at low 
Ay in the WD01-7 model. This is a result of the higher gas tem- 
perature due to the stronger PE heating efficiencies of that dust 
model. 

The right panel in Fig.[2T|shows the isolated effect of the dif- 
ferent H2 formation treatment. All models in the right panel were 
computed using the WD01-7 dust properties regarding the radia- 
tive transfer and the PE heating. The C + -CO transition is similar 
across all shown models, since the FUV intensity, which con- 



trolls the photo-dissociation of CO, is the same across all mod- 
els. However, the CO abundance at lower Ay shows very large 
differences. The models using the CT dust formation scheme 
produce more CO at the cloud edge compared to the SD89 H2 
formation. The CT models provide a stronger H2 formation heat- 
ing contribution. The additional heating term produces a signifi- 
cant population of hot CO in the outer layers of the model clump. 

The production of hot, i.e., strongly excited, CO is visible 
in the spectral line emission of the model clumps. In Fig. [22] 
we show the corresponding total CO line emission of the model 
clumps shown in the right panel in Fig. [21] For transitions higher 
than J=8-7, the line-integrated intensities start to show signif- 
icant differences for the three models. In this particular case, 
the strength of the high-/ CO transitions is proportional to the 
total H2-forming grain surface. The outer layers of the model 
clumps contribute significantly to the total CO line emission 
of this model. As a consequence, any modeling with a large 
surface of the clouds, e.g., assuming a clumpy structure, will 
show much stronger high-/ CO emission lines compared to non- 
clumpy model attempts. 

As an additional example of chemical species that predomi- 
nantly form in the outer regions of molecular clouds and thus are 
affected by the H2 formation treatment we show in Fig. [23] the 
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Fig. 23. CH and CH + abundances (left) and column densities (right) for a model clump of M = 10 M , n — 10 5 cm~ 3 , and 
X = 1000. The different lines show the effect of different H2 formation treatment. All models were computed using the WDO 1-7 dust 
properties regarding the radiative transfer and the PE heating. 
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Fig. 24. OH and OH + abundances (left) and column densities (right) for a model clump of M = 10 M , n — 10 5 cm~ 3 , and 
X = 1000. The different lines show the effect of different H2 formation treatment. All models were computed using the WDO 1-7 dust 
properties regarding the radiative transfer and the PE heating. 



CH and CH + abundances and column densities and in Fig. [24] 
the OH and OH + abundances and column densities. All these 
species are significantly affected by the strong H2 formation 
heating effect in the outer parts of the model clump. The col- 
umn density effect is weakest for the CH because of the abun- 
dance peak at Ay « 1. The other three species show a stronger 
effect on the column density because of the weak formation in 
deeper clump regions. Recent Herschel observ ations show high 
column densities of these surface tracers (e.g. lOin et alj|2010l : 
iFalgarone et al.l20Toh . which were difficult to reproduce in exist- 
ing PDR models. A more efficient H2 formation modeling with 
the accompanying formation heating effect may help to resolve 
the discrepancy between models and observations. 



4.2.1. Full model SED 

In Fig. [25] we show the full continuum and spectral line flux 
(in units of Jy) emitted by a model clump of M = 10 M Q , 
n = 10 5 ctrT 3 , and^f = 1000. The two panels show the results 
for different dust models. The top panel shows the result for the 
dust model with the largest effective surface contributing to the 
gas heating, the WD01-7 model, equivalent to a Ry = 3.1. The 
strong gas heating is visible in the large population of intense 
high-/ CO emission lines with J> 10 stemming from the outer 
layers of the cloud (see Table[5]).Q In the bottom panel we show 
the result for a MRN dust model, visible through the lack of PAH 
emission features in the NIR. 



7 The level population of the CO energy levels was computed up to 
J=50. 
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Table 5. Line emission J Idv for the three models shown in Fig. [25] 





line 


A 


WD01-7 
w/o PAHs 


WD01-7 
with PAHS 


WD01-25 
w/o PAHS 


WD01-25 
with PAHS 


MRN 






[}im\ 


[erg/s/cm 2 /.sr] 


[erg/s/cm 2 /jT] 


[erg/s/cm 2 /.sr] 


[erg/s/cm 2 /.sr] 


[erg/s/cm 2 /.sr] 


^TIR 






0.290 


0.290 


0.286 


0.280 


0.278 


[Cn] 


2 Pl/2 - 2 Pl/2 


158 


7.5 x 10~ 4 


7.5 x 10- 4 


9.2 x 10~ 4 


9.3 x 10- 4 


7.7 x 10~ 4 


[Oi] 


3 'Pi- 3 Pi 


63 


6.8 x 10" 3 


6.6 x 10" 3 


5.7 x 10 -3 


6.0 x 10~ 3 


5.5 x 10~ 3 


[Oi] 


3 P - 3 Pi 


146 


5.0 x 10~ 4 


5.1 x 10- 4 


4.2 x 10~ 4 


4.9 x 10~ 4 


4.1 x 10- 4 


[Ci] 


3 Pi - 3 Po 


609 


1.2 x 10~ 6 


1.2 x 10~ 6 


1.4 x 10~ 6 


1.3 x 10~ 6 


8.2 x 10~ 7 


[Ci] 


3 Pl~ 3 P l 


370 


3.9 x 10~ 6 


3.9 x 10- s 


4.6 x 10~ 6 


4.2 x 10~ 6 


2.4 x 10- 6 


CO 


(1-0) 


2601 


1.3 x 10~ 7 


1.3 x 10~ 7 


1.2 x 10- 7 


1.3 x 10~ 7 


1.5 x 10- 7 


CO 


(4-3) 


650 


7.3 x 10~ 6 


7.4 x 10- s 


7.0 x 10~ 6 


7.5 x 10~ 6 


8.9 x 10- 6 


CO 


(7-6) 


372 


1.1 x 10" 5 


1.1 x 10" 5 


1.3 x 10" 5 


1.5 x 10~ 5 


1.9 x 10" 5 


CO 


(10-9) 


260 


2.7 x 10~ 7 


5.0 x 10- 7 


2.8 x 10~ 7 


4.5 x 10~ 7 


7.5 x 10- 7 


CO 


(15-14) 


174 


5.6 x 10- 8 


1.6 x 10~ 7 


2.4 x 10"'° 


2.5 x 10~ 8 


8.4 x 10-" 


CO 


(20-19) 


130 


1.2 x 10~ 8 


3.7 x 10~ 8 


8.7 x 10 12 


7.3 x 10-" 


1.7 x 10-" 



100 
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Fig. 22. Total line integrated clump-averaged intensities of 12 CO 
transitions for a model clump of M — 10 M Q , n — 10 5 cm -3 , 
and^ = 1000. The different lines show the effect of different 
H2 formation treatment. All models were computed using the 
WDO 1-7 dust properties regarding the radiative transfer and the 
PE heating. 

In Table [5] we summarize the line emission of the WD01- 
7 and WDO 1-25 (with and without H2 formation on PAHs), as 
well as the MRN models quantitatively. The differences between 
the WDO 1-7 and WDO 1-25 models are limited to the line emis- 
sion, the continuum emission is only marginally different. The 
WDO 1-25 models show weaker high-/ CO emission because of 
the smaller grain surface and correspondingly a weaker H2 for- 
mation heating. The high-/ CO emission of the MRN model falls 
between the WDO 1-7 and the WDO 1-25 dust models without H 2 
formation on PAHs. Allowing for H2 formation on PAHs will 
add a strong heating contribution from the formation processes 
on the large PAH surface, exciting the high-/ CO transitions by 
a large factor. 



5. Summary 

We revised the treatment of the dust in the KOSMA-t PDR 
model code to achieve a consistent description of the dust re- 
lated physics in the code. The computation of dust properties 
in numerical models of photo-dissociation regions afflicts the 
chemical and physical structure of a model cloud via multiple 
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Fig. 25. Full continuum and spectral line flux emitted by a model 
clump of M = 10 M , n = 10 s cm 4 , and x = 1000. The two 
panels show the results for different dust models. Top: WD01-7 
dust model, equivalent to a R y = 3.1. This is the dust model with 
the largest grain surface contributing to the heating. Bottom: 
MRN dust mode. This is a dust model with a slightly smaller 
grain surface com pared to WDO 1-7 . The P E heating was calcu- 
lated according to lBakes & Tielensl (11994 . 



effects. The major areas where dust properties play an important 
role are 
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1. the optical dust properties, i.e., their influence on the radia- 
tive transfer 

2. the dust temperature (neglecting non-equilibrium heating for 
the very small grains), which influences surface chemistry, 
most importantly the formation of molecular hydrogen 

3. the heating and cooling capabilities of dust grains via photo- 
electric heating and gas-grain collisional cooling. 

The effect of changes in these areas on the physical and chemi- 
cal structure of a model can be profound and we described their 
respective influence in detail. 

We notice two opposite effects: Increasing the fraction of 
small grains steepens the slope of the extinction curve so that 
the FUV intensity drops faster when going into the cloud, which 
leads to a temperature decrease at high optical depths. However, 
it also provides more surface for PE heating and H2 formation, 
thus increasing the heating efficiency. As a consequence, the gas 
becomes hotter in the outer layers, but colder inside the cloud 
when the fraction of small grains and PAHs is increased. Only 
tracers of the outer layers are able to reliably measure the gas 
heating efficiency. In contrast to common wisdom, this is not 
even true for ionized carbon [CII], which traces both regimes. 

The most influential modification of the code is the treat- 
ment of H2 formation on the surface of dust grains. Allowing for 
chemisorption, i.e., including the Eley-Rideal effect into the for- 
mulation of the Information efficiency, significantly increases 
the overall formation rate of molecular hydrogen. If a significant 
fraction of the dust consists of very small grains and PAHs, their 
very large contribution to the total grain surface leads to an addi- 
tional enhancement of the total H2 formation. If the molecule 
formation is enhanced, the formation heating is enhanced as 
well, providing a profound heating contribution in the outer lay- 
ers of the model clumps, which significantly affects the local 
chemistry. 
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Appendix A: Isotropic vs. non-isotropic scattering 

The MCDRT code assumes isotropic scattering, i.e., a mean 
scattering angle g = (cos©) = when calculating the 
frequency-dependent FUV radiative transfer. However, taking 
into account the detailed material properties and size distribution 
of interstellar dust results in different values of g. Reasonable 
values for int erstellar dust in the FUV range fall around g x 0.7 
and w « 0.4 ihi & Drainell2001h . In the following we show that 
even in this case of strong forward-scattering the isotropic scat- 
tering treatment poses a clear improvement compared to the pure 
forward-scattering cased used in our previous model. 

To quantify the error we used the s pherical harmonics 
method, presented bv lFlannerv et al.1 (1 19801) to calculate the pen- 
etration of UV radiation of a spherical cloud. Applying their 
method, we calculated the attenuation of the mean intensity as 
a function of optical depth for a given asymmetry fac tor g and 
albedo at. For more details on the method see also Robergel 
(l!983l) : lLe Petit et al.l J2006l) : lGoicoechea & Le Bourlotl(l2007l) . 

In Fig. IA.1I left panel, we plot the solution of the mean in- 
tensity for a) = 0.4, g = 0.7 normalized to the incident radiation 
field for a homogeneous spherical cloud with an optical depth to 
the center t c = 20 (black curve). The dashed blue curve shows 
the solution for the same albedo, but for isotropic scattering. 
This represents the solution of the MCDRT code. The dotted red 
curve gives the result for pure forward-scattering, i.e., the case 
exp(-(l - w)t). This is equivalent to the previous approximation 
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Fig. A.l. Left panel: Mean intensity versus optical depth r normalized to the incident specific intensity Iq. The solid black curve 
shows the depth dependence of the mean intensity assuming dust parameters a> = 0.4, g = 0.7, approximately describing the 
WD01-7 dust sort. The dashed and dotted curves denote mean intensities using the same albedo to, but assuming isotropic and 
pure forward-scattering, respectively. Right panel: Absolute error | log(y,//o) - \og(Jklh))\ of assuming isotropic g = and pure 
forward-scattering g = 1 with respect to more realistic values <d = 0.4, g = 0.7 (dashed and dotted lines, respectively). 



in KOSMA-t. The right panel shows the corresponding error 
when assuming isotropic or pure forward-scattering. It is appar- 
ent that the assumption of isotropic scattering is an improvement 
over our previous approach up to optical depths of about 16, only 
for very high t the forward-scattering is the solution closer to the 
g = 0.7 case. However, the FUV influence at t > 10 is negligi- 
ble. Presently, we therefore sacrifice the dark cloud accuracy in 
favor of a more accurate description of the cloud surface, but in 
a future work we plan to update the MCDRT code to properly 
account for non-isotropic scattering. 



Appendix B: Rescaling of photo-rates to account 
for different dust properties 

In the following we discuss the uncertainties of the proposed 
rescaling jj — > for different dust models D. 

One important factor is the spectral shape of the assumed 
FUV radiation field. In Fig. IB.ll we plot the photo-dissociation 
cross sections cr of two different molecules (CH, in the left 
panel, SiH + in the right panel) together with the FUV radiation 
field (thick black line: unshielded Draine field, colored lines: 
FUV field inside the clump at 95% of the total radius for the 
dust models WD01-7, 21, and 25). The total photo-dissociation 
cross section is calculated using Eq. [13] The different wave- 
length behavior of cr in the two panels demonstrates why the 
photo-dissociation rates of different species are attenuated dif- 
ferently when going deeper into the cloud. The right panel also 
emphasizes that the choice of the spectral shape of the illumi- 
nating FUV radiation field may influence the photo-dissociation 
rate of some species significantly; the choice between the two 
FUV fields results in an unshielded photo-dissociation rate dif- 
ference of a factor two. However, differences in the attenuation 
behavior between different dust models will still be dominated 
by wavelengths shorter than 2000 A because A(A)/A(V) varies 
only weakly for A > 3000 A. 




5 10 15 

A v 



Fig. B.2. Comparison of the explicitly integrated photo- 
dissociation rates of CH 2 for the WD01-7, WD01-21, and 
WD01-25 dust models (circles, squares, triangles, respectively) 
and the corresponding fits to the rates according to Eq. [14] The 
different line styles denote fits up to different maximum visual 
extinctions (solid: Ay.max = 3, dashed: Ay, m ax = 10, dotted: 
^v.max = 15, and dot-dashed: A v , ma x = 20). 



Another uncertainty results from the choice of Ay.max up to 
which the jjju from Eq.[14]is fitted to the values Tj(Ay)/Tj(0), 
In Fig. IB. 21 we show the Tqh (^v)/Fch,(0) for three dust mod- 
els: WD01-7, WD01-21, and WD01-25. The lines shows the 
corresponding exp(-y^ )( Av) fits for different Ay >m ax- The range 
of yf fu is 1.92-2.24, 1.48-1.71, and 1.14-1.3 for the three dust 
models, equivalent to a 10-15% uncertainty. We calculated fits 
to yffn for all available species and for all dust models to 

^V.max = 10. 

The third source of uncertainties is the functional depen- 
dence used in the fit through Eq.[l6] We performed least-squares 
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Fig. B.l. Wavelength dependence of the assumed UV illumination and the photo-dissociation cross sections. The thick black line 
shows the unshielded empty-space Draine FUV field, the dashed line shows an extension of the Draine field longward of 2000 
A suggested by Ivan Dishoeck & Blackl dl982l) . The blue, green and red lines show the FUV radiation field inside the reference 
cloud at 95% R tot for the WD01-7, WD01-21, and WD01-25 dust models, respectively. Left panel: The blue filled curve shows the 
continuum absorption cross section cr cont of CH 2 that leads to dissociation, the black dots with drop lines show the corresponding 
line absorption cross section <xij ne . Right panel: Same as left panel but for the SiH + photo-dissociation cross sections. 



fits to polynomials y® = p{jj, l/^^) U P to order 4 to de- 
rive a suitable scaling relation between the input parameters jj, 
describing the wavelength-dependent attenuation of the photo- 
dissociation rate of species j with Ay, l/^ 6() ^ describing the 
wavelength-dependent attenuation of the FUV radiation field 
depending on the dust model D, and the output parameter y 1 ?. 
Figure lB~3l visualizes Eq. [16] (black grid) and the explicitly cal- 
culated fits yfj it (colored dots) for the 25 WD01 dust models 
for all species with available cross sections. The data points are 
colored according to their deviation from Eq.[T6l 

Seven species, CH + , SH + , OH + , HCO + , CO, Ot and SiO, 
deviate significantly from our Eq. [16] with residual >0.4. By 
shifting their yj to a corrected fj we are able to compensate for 
these deviations. Figure IB .41 shows the fitted y?, fl of the three 
strongest outliers and the corresponding y® from Eq. [16] using 
the original y y - (dashed lines) and the corrected yj (solid lines). 

The deviations \jjju — yf\ of the seven outlying species 
are large compared to the rest of the 63 species that we calcu- 
lated (compare Fig. IB. 3b . In the following we will discuss the 
strongest outliers individually^. 



Draine FUV field at longer wavelengths should not affect 
the final photo-dissociation rate. However, computing the un- 
shielded photo-dissociation rate coefficient a for a Draine field 
yields 1.27 x 10'" s' 1 comparable to 1.1 x 10~ 11 s -1 from 
Ivan DishoecketaLl (12006]) . while UDfA lists 1 x 10 12 s 1 . The 
reason for the mismatch is unclear, al l sources refer to the same 
photo-dissociation cross section by I Saxon &Liul (Il986h . We 
also find inconsistent values of y published, between 1.8 and 
3.5. We find a best fit to Eq.[[6]with yf = 2.76. 



CH + 



For CH + dev i ations from Eq. [16] of up to 0.6 were found. 
IWoodall et all d2007l) give a = 2.5 x 10' 1Q s' 1 and y = 2 .5 re- 
ferring to lRoberge et al.l d 19911) . Ivan Dishoeck et al ] d2006l) give 
a = 3.3 x 10' 1Q 
iKirbv et al.l (fl980l) . We find y. f . 

Am these lower values of y lead to photo-dissociation rates 3-6 
orders of magnitude stronger compared to rates calculated with 
the higher y from above. The reason for the y offsets is unclear. 



1 and y = 2.94 using cross sections from 
D - 1.5 -2.4 and yf = 2.11. At 



HCO + 

HCO + showed the strongest maximum deviation of 1.12 from 
th e plane defined by Eq. [T6l UDfA gives y = 2.0 referring 
to Ivan Dishoeck et aD ([2()()6) who give a value of y = 3.32 
in their Tab. 2. Unfortunately, all photo-dissociation rate coeffi- 
cients from UDfA that refer to van Dishoec k et al. I (120061) are in- 
consistent with their published numbers and should be replaced 
with more recent numbers. Their paper was still in preparation 
when IWoodall et ail d2007l) was published, so the mismatch is 
possibly due to pre-publication updates of the final values. 

OH + 

OH + showed the second-strongest deviations of up to 0.8 from 
Eq. [16] It possesses a strong absorption continuum only for 
wavelengths shorter than 1100 A so that extensions to the 

8 We ignore deviations for CO b ecause photo-dissociation of CO is 
treated separately in any PDR code dvan Dishoeck & Blacklll988h 



SH + 

SH + has comparable val ues of a, and y in UDfA and in 
van Dishoeck et a D d2006l) . Small differences in a can be at- 
tributed to a strong absorption line at 3100 A, which leads to 
a stronger unshielded photo-dissociation rate w hen applying a 
10000 K black-body radiation field as used by Woodall et al. 
d2007h . We find y° = 1.2-1.8 and yf = 1.39. Differences 
of y from the literature and our calculations most likely re- 
sult from different assumed dust properties and fits to lower 
maximumA V; max- 

SiO and O* 

Both molecules show maximum deviations of 0.4 from 
Eq. [161 The unshielded S iO rate coefficients from UDfA and 
Ivan Dishoeck et ail d2006l) differ by a factor 16 (see comment on 
HCO + ) while their yj are the same. 
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6 

Fig. B.3. Visualization of Eq. [16] (black grid) and the explicitly calculated fits y?„ i( (colored dots) for the 25 WD01 dust models 
for all species with available cross sections. The data points are colored according to their deviation from Eq. [16] Green and blue 
denotes deviations of less than 0.05 and 0.1. Red points signal deviations exceeding 0.1. All deviations are smaller than 0.22. 



In Fig. IB. 51 we show a histogram of the residuals yfj it - 
yR after removing the outlying species. The standard de- 
viation is cr — 0.05, the maximum absolute residual 
is 0.22. Photo-dissociation cross sections are available for 
a number of species that are n ot yet included in UDfA 
(van Hemert & van Dishoeck 2008). In Tab. IB.ll we give a and y 
from http://home.strw.leidenuniv.nl/~ewine/photoas well as our 

7j- 



Appendix C: H 2 formation rate 

We followed the model of H2 formation on interstel- 
lar dust grains via physisorption and chemisorpti on from 
Cazaux&Tielens dCazaux & Tielensll2002ll2004l2010l) . 

For a given dust type, the total H2 formation rate is given by 
equation [19] (omitting the subscript index ;' in the following) 



Rh 2 = -nHV H nrfcr rf eH 2 5 H , 
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Fig. B.4. Comparison of the y?^.. ( for the three strongest outliers 
CH + , OH + , and HCO + and the corresponding scaling behavior 
of Eq. [16] using the original jj (dashed line) and the corrected 

where the formation efficiency is given by equation l22l 
e H2 = (* + 1 + HY X $ 

+ l + (C1) 
\2j8h 2 a Pc ) 
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Fig. B.5. Histogram of the residuals y?j it - yf after removing 
outliers. 

Table B.l. Rate coefficients and corrected jj for species not in- 
cluded in UDfA. PI and PD denotes photo-ionization and photo- 
dissociation. 



process 


species 


a(s-') 




7j 


Ji 


PI 


CH 4 


6.8 x 10 


12 


3.94 


4.13 


PD 


CS 2 


6.1 x 10 


9 


2.06 


1.75 


PI 


CS 2 


1.7 x 10 


9 


3.16 


2.68 


PD 


H0 2 


6.7 x 10 


10 


2.12 


1.81 


PD 


N,0 


1.9 x 10 


9 


2.44 


2.02 


PI 


N 2 


1.7 x 10 


10 


3.93 


4.06 


PD 


NH + 


5.4 x 10 


11 


1.64 


1.32 


PI 


N0 2 


1.5 x 10 


HI 


3.33 


2.86 


PD 


o 3 


1.9 x 10 


9 


1.85 


1.48 



We set (jj.F)/(2f3u^) to zero to make sure that newly formed H2 
molecules are a ble to leave very cold dust surfaces, equivalent to 
equation (13) in lCazaux & Tielensl (120021) . We can approximate 
£ and PH p /a PC by 



f =1 + 



vh c exp 



1.5E 



He 



1 + 



'He 



K2F)' 




(C.2) 
(C3) 



-1 



8 \JnT d exp i-2a Es l \ exp y£ Hc - E Hp 



Equation IC.3I takes i nto account the corrections from 
ICazaux & Tielensl d2010l) . Table dC.U gives the applied model 
parameters for silicate and graphite surfaces. 

Appendix D: Photoelectric heating rates 

IWeingartner & Draind d2001al) showed that the photoelectric 
heating rate for the WD01 grain size distributions is fairly well 
reproduced by the parametrization 

Cq + CiT^ 



r pc = Gn H 



1 +C 2 ^ C5 [1 +C 3 iA Cf >] 



10" 26 ergs" 1 cm" 3 , 



Table C.l. Model parameter for silicate and graphite surfaces. 
For more detail s about the determination and calculation of these 
parameters, see lCazaux & Tielensl ( l2002f) 



Parameter 


silicate 


graphite 


E H , [Kf 


320 


520 


n" 


0.3 


0.4 


Es [K] c 


110 


260 




450 


800 


En c [K]' 


30000 


30000 


j2m(E M -E s ) f 

"V £ f 


14.4 


14.0 


VH 2 [S-'] S 


3x 10 12 


3x 10 12 


v Hc [s-']" 


1.3 x 10 13 


1.3 x 10 13 


F [mLys -1 ]' 


10 


10 



Notes. (a) Desorption energy of H 2 . (A) Fraction of newly formed 
H 2 that stays on the surface. (<) Energy of the saddle point between 
a chemisorbed and a physisorbed site. {d) Desorption energy of ph- 
ysisorbed H. w Desorption energy of chemisorbed H. (fl Product of 
the width and height of the energy barr ier between a phys i sorbed and a 
chemisorbed site. For more details see lCazaux & Tielensl d2004l) (g) vi- 
brational frequency of H 2 in surface sites. {h) vibrational frequency of H 
in surface sites. (,) Flux of H atoms in monolayers per second. 




Fig. D.l. Fit results to the photoelectric heating rate, calculated 
for an extended range of G y/T/n e . The different lines show re- 
sults for gas temperatures 100 K, 1000 K, and 10000 K. 
6 ergs~'. 



given in units of 10 21 



r pc is 



with if/ — (G y/T)/n e , where T is the gas temperature in Kelvins 
and 4> is in units of K'^cm 3 , and G = l-Hx- Numerical val- 
ues for the parameters C, in E q. dD.U are given in Table 2 in 
IWeingartner & Draind d2001al) for each dust distribution. The 
rate of cooling due to charged particle collisions can be approx- 
imated by 



A Rr = n e tinT 



Do+Di/x 



exp D 2 +D3X- D4X" 



(D.2) 



xlO 28 erg s 1 cm 3 , 



where x = In (G Vr jn e ). Values for the A for all WD01 
grain di stributions are given in Table 3 in lWeingartner & Draind 
(I2001al) . 

Equations dD.U and (lD.2t are fairly accurate when 10 3 K < 
T < 10 4 K and 10 2 K 1/2 cm 3 < ifr < 10 6 K 1/2 cm 3 . However, the 
parameter range in PDR models allows i[s to rea ch values as 
high a s 10 9 K^ 2 cm 3 and as low as 10~ 2 K^ 2 cm 3 . IWeingartnerl 
(120091) provided updated calculations of r pc and A gr for the 
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Table D.l. Photoelectric heating parameters (equation [25). See IWeingartner & Drainel (1200 1 ah for detailed descriptions of the dust 
size distributions and applied physical parameters. 



Rv 


b c 


Distr. 


Rad. Field 


Co 


Cj 


c 2 


c 3 


c 4 


c 5 


c 6 


c 7 


3.1 


6.0 


A 


B0 


8.36780 


0.20755 


0.00969 


0.01354 


0.41047 


0.62072 


0.53099 


0.00589 


4.0 


4.0 


B 


B0 


5.67140 


0.13563 


0.00994 


0.01568 


0.41122 


0.60659 


0.53265 


0.00382 


5.5 


3.0 


B 


B0 


5.67140 


0.13563 


0.00994 


0.01568 


0.41122 


0.60659 


0.53265 


0.00382 



Table D.2. Collisional cooling parameters (equation [26). See IWeingartner & Drainel (1200 lal) for detailed descriptions of the dust 
size distributions and applied physical parameters. 



Rv 


b c 


Distr. 


Rad. Field 


Do 


Dj 


D 2 


£>3 


D 4 


D 5 


D 6 


Z> 7 


Dg 


Dg 


Dio 


£>„ 


3.1 


6.0 


A 


BO 


0.493 


0.556 


0.031 


0.842 


287 


0.360 


3.640 


-0.098 


0.642 


3.25 x IO" 6 


2.080 


-1.500 


4.0 


4.0 


B 


B0 


0.493 


0.362 


0.038 


0.803 


287 


0.360 


2.460 


-0.035 


0. 


2.35 x 10~ 5 


1.870 


-1 


5.5 


3.0 


B 


B0 


0.484 


0.720 


0.118 


0.597 


100 


0.710 


1.240 


-0.030 


0.603 


1.11 x 10~ 6 


2.150 


-2.720 



lo 4 F 2 " 




i<r 2 io-' 10" 10 1 10 2 10 2 10 4 10 s 10 6 10 7 10 s 10 4 

GVT /n, IK" 2 ™ 2 ] 



Fig. D.2. Fit result to the net cooling rate due to collisions with 
charged particles, calculated for an extended range of G Vr / 
The different lines show results for gas temperatures 100 K, 
1000 K, and 10000 K. F pc is given in units of 10 -26 ergs -1 . 

dust distributions WD01-7, WD01-21, and WD01-25 up to if/ = 
lO'K^cm 3 . 

We performed numerical fits to the updated photoelec- 
tric heating rates and to the cooling rate due to collisions 
with charged particles according to the parametrizations 1251 
an d l26l The heating a nd cooling rates have been calculated 
by IWeingartnerl (120091) with an extended parameter range of 
10" 2 K 1/2 cm 3 < (/> < 10 9 K 1/2 cm 3 . The details of t he computa- 
tion are described in IWeingartner & Drainel d2001ah . To achieve 
a good fit over the full extended parameter range it was neces- 
sary to modify the alg ebraic form of the parametrizatio ns relative 
to its original form in IWeingartner & Drainel (l2001al) . The final 
form is given in equations 1251 and l26l In Figs. ID.ll and lD.2l we 
show the fits to calculations of F pc and A gl for gas temperatures 
of 100 K, 1000 K and 10000 K to demonstrate the quality of the 
parametrization. The C, and D, parameters are given in Table lD.ll 
and lD.2l and will provide a fairly accurate reproduction of the to- 
tal photoelectric heating rate F tot between 10 K < T < 10000/T 
and 10- 2 K 1/2 cm 3 < if/ < 10 9 K 1/2 cm 3 . 
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